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Abstract 

Based on the concept of self-decomposable random variables we discuss the 
application of a model for a pair of dependent Poisson processes for energy 
facilities. Due to the resulting structure of the jump events, we can see the 
self-decomposability as a form of cointegration among jumps. In the con¬ 
text of energy facilities, the application of our approach to model power or 
gas dynamics and to evaluate transportation assets seen as spread options is 
straightforward. We study the applicability of our methodology first assum¬ 
ing a Merton market model with two underlying assets; in a second step we 
consider price dynamics driven by an exponential mean-reverting Geometric 
Ornstein-Uhlenbeck plus compound Poisson that are commonly used in the 
energy field. In this specific case we adopt a price spot dynamics for each 
underlying that has the advantage of being treatable to find no-arbitrage con¬ 
ditions. In particular we can find close-form formulas for vanilla options so 
that the price and the Greeks of spread options can be calculated in close 
form using the Margrabe formula [6] (if the strike is zero) or some well known 
approximations as in Deng et al. [10] . 

*The technique here discussed does not reflect UGC view. 

Icufaro@ba.infn.it 

I piergiacomo .sabino@uniper. energy 


1 



N Cufaro Petroni and P Sabino: Cointegrating jumps for energy facilities 


2 


1 Introduction and Motivation 

Several research studies have shown that the spot dynamics of commodity prices 
is subjected to mean reversion, seasonality and jumps (see for example Cartea and 
Figueroa m)- In addition, some methodologies have also been proposed to take 
dependency into account based on correlation and cointegration. However, these 
approaches can become too mathematically complex or non-treatable when leaving 
the Gaussian-Itd world. 

In this paper we address the problem of dependency in the 2-dimensional case 
and start considering 2-dimensional jump diffusion processes with a 2-dimensional 
compound Poisson component. We then introduce an intuitive approach to model 
the dependency of 2-dimensional Poisson processes based on the self-decomposability 
(see Cufaro Petroni |8], Cufaro Petroni and Sabino [9], Sato [12]) of the exponential 
random variables used for its construction. In Cufaro Petroni and Sabino [S] we have 
seen that given two independent exponential rv’s H, Z ~ S(A), and a 0-1 Bernoulli 
rv B{1) ~ iB(l, 1 — a) with a = P {5(1) = 0}, then also the rv defined as 

X = aY + Za Za = B{l)Z (1) 

is an exponential C(A) resulting in a weighed sum of Y and Z where 0 < a < 1 is the 
deterministic weight of Y, while 5(1) is the random weight of Z. In other words X 
is nothing else than the exponential Y down a-rescaled, plus another independent, 
but intermittent with frequency 1 — a, exponential Z. It is apparent on the other 
hand that, by construction, X and Y are not independent and it is possible to show 
that a also represents precisely their correlation coefficient. This result is a direct 
consequence of the self-decomposability of the exponential laws. 

As a matter of fact we could also produce pairs of a-correlated exponentials 
X ~ S(A) and Y' ~ with different parameters by reformulating the previous 
relation as 

A = 7F' + Z, (2) 

where 7 = pa/A, so that 0 < 7 < p Considering then X and Y' as two random times 
with a positive random delay Za, the mathematical concept of self-decomposability 
can help describing their co-movement and can answer some common questions 
arising in the hnancial context: 

• Once a hnancial institution defaults how long should one wait for a dependent 
institution to default too? 

• A market receives a news interpreted as a shock: how long should one wait to 
see the propagation of that shock onto a dependent market? 

• If different companies are interlinked, what is the impact on insurance risk? 

Questions like the ones above are covered by the special case 7 > 1. Our model 
is then rich enough to describe cases where the second random time event occurs 


N Cufaro Petroni and P Sabino: Cointegrating jumps for energy facilities 


3 


before the first one. Similar results based on linear structure of exponential rv can 
be found in Iyer et al. [15] whose purpose was to model a multi-component reliability 
system. 

It is worthwhile to notice that this bivariate exponential model implies a copula 
function (see Cufaro Petroni and Sabino [9]) that is neither chosen upfront nor 
whose parameters are estimated from market data: here the copula function does 
not dehne the model, but rather the opposite. 

Based on the self-decomposability of the exponentional rv’s we are able to 
construct 2-dimensional Poisson processes with dependent marginals (see Cufaro 
Petroni and Sabino m)- Because of the relationship among random times, the two 
Poisson processes can be seen linked with a form of cointegration between their 
jumps. 

In the context of energy facilities, the application of our approach to model price 
dynamics and to evaluate transportation assets seen as spread options is straight¬ 
forward. To this purpose, we consider the German EEX and French Powernext 
power markets and assume that each spot price dynamics is driven by an exponen¬ 
tial mean-reverting Geometric Ornstein-Uhlenbeck (GOU) plus compound Poisson. 
In this specihc case we adopt a stochastic dynamics of the spot prices that is slightly 
different from the one in Gartea and Figueroa [T] with the advantage of being more 
treatable to hnd no-arbitrage conditions. In particular we can End close-form for¬ 
mulas for vanilla options hence the price and the Greeks of spread options can be 
calculated in close form using the Margrabe formula [6] (if the strike is zero) or some 
well known approximations as in Deng et al. [101, Pellegrino and Sabino ini and Pel¬ 
legrino [7] . In any case our approach implies an explicit algorithm for the simulation 
of the dependent Poisson processes and can be used in Monte Garlo simulations. 

Finally, we compare the results obtained by our approach to the ones obtained 
assuming that the two compound Poisson processes are independent or contain a 
common Poisson component. 

The extension to the multi-dimensional case will be the goal of future studies 
as well as the extension to different dynamics other than Poisson. However, under 
the assumption that only two underlyings have jump component, the price and the 
Greeks of spread options can be obtained by the moment-matching methodology 
proposed in Pellegrino and Sabino [T6] . 

The paper is organized as follows. Section [2] summarizes the results for bi- 
dimensional Poisson process introduced in Gufaro Petroni and Sabino [9]. In section 
Owe consider bi-dimensional jump diffusion processes having a Geometric Brownian 
Motion (GBM) and GOU diffusive component. We also apply our methodology to 
the 2-factor Schwartz-Smith model [T3| with jump diffusion where we hnd analytical 
solutions for vanilla options as well. Section 0] presents the risk neutral formulas for 
plain vanilla and spread options given the price dynamics introduced in Section[3]and 
given the different types of bi-dimensional Poisson components. Section 0] illustrates 
our approach with practical examples: we hrst assume a pure GBM plus jump model 
(Merton model) and compare the results obtained by our approach to the ones 
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obtained assuming that the two compound Poisson processes are independent or 
contain a common Poisson component. In a second step we estimate the parameters 
of a 2-dimensional GOU plus jumps dynamics to model the EEX and Powernext 
day-head prices. Finally we compare the price of an interconnection between these 
two locations assuming the three types of Poisson conhguration mentioned above. 
Section [6] concludes the paper with an overview of future studies and possible further 
applications. 

2 Dependent Poisson processes 

A law with density (pdf) f{x) and characteristic function {chf) (p{u) is said to be 
self-decomposable (sd) (see Sato |T2] and Cufaro Petroni [8]) if for every 0 < a < 1 
we can End another law with pdf ga{x) and chf Xa{.u) such that 

(p{u) = (p{au)xa{u) 

This dehnition selects an important family of laws with many relevant properties. 
Remark however that, while a sd Xa{u) can be explicitly expressed in terms of (p{u), 
its corresponding pdf ga{x) can not be given in a general, elementary form from 
f{x). We will also say that a random variable {rv) X is sd when its law is sd: 
looking at the dehnition this means that for every 0 < a < 1 we can always End two 
independent rv’s Y (with the same law of X), and Za with pdf ga{x) and chf Xa{u) 
such that in distribution 

X = aY + Za 

We can look at this, however, also from a difierent point of view: to the extent that 
for 0 < a < 1 the law of Za is known, we can deEne the rv 

X = aY + Za 

which by self-decomposability will now have the same law of Y. It would be easy 
to show that a also plays the role of the correlation coefEcient between X and Y, 
namely 

rxY = a 

It is well known, in particular, that the exponential laws Ci(A) with pdf and chf 
fi{x) = Ae"^^l„>o (pi{u) = ^ . 

A — tU 

are a typical example of sd laws (see Sato [T2] I. 

It is possible to show now (Cufaro Petroni and Sabino [9]) that the law of Za is 
a mixture of a law <5o degenerate in 0, and an exponential Ci(A), namely 


Za ~ -|- (1 — a)Ci(A) 
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and this entails that Za can be taken as the product of two id rt>’s: Z ~ Si(A), and 
B{1) 03(1,1 — a) (a Bernoulli with a = P {-B(l) = 0}), which are also independent 

from V, namely 

Za = B(l)Z 

In summary, given two exponential rv^sY ~ Si(A) and Z ~ C:i(A), and a Bernoulli 
B{1) rsj 03(1,1 — a) (all mutually independent) the rv 

X = aY + B{l)Z 

is again an exponential Ci(A) dehned as the weighed sum of Y and Z: while a is 
the deterministic weight of Y, the weight of Z is random and is represented by 
another (independent from both Y and Z) 0-1 Bernoulli rv B{1) ~ 53(1,1 — cl)- 
In other words X is nothing else than the exponential Y down a-rescaled, plus 
another independent, but intermittent with frequency 1 —a, exponential Z. The self- 
decomposability of the exponential laws ensures then that, if both the parameters 
of Y and Z are A, also X marginally is an Si (A) for every 0 < a < 1. It is apparent 
on the other hand that, by construction, X and Y are not independent and it is 
easy to show that a represents their correlation coefficient. 

Remark that if Y' ~ Si(/i), then aY' ~ Si (^) for every a > 0, and hence in 
particular 

Y = ^Y' ^ Si(A) 

A 

As a consequence we could also state the self-decomposability by means of exponen¬ 
tial m’s with different parameters X ~ Si (A) and Y' ~ Si(/i) because of course we 
have 

X = aY + Za = ^Y' + Za 

provided that 0 < a < 1 . 

As initially suggested in Iyer et al. [15], we take now a sequence of iid rv's 
Xk = aYk + Rfc(l)Zfc k = 1 , 2 ,... 

in such a way that for every k\ Xk,Yk,Zk are Si(A), 5^(1) is 53(1,1 — a), and 
Yk, Zk, Bk{l) are mutually independent. Add moreover Xq = Yq = Zq = 0, P-a.s. 
to the list, and then dehne the point processes 

n 

Tn = 5^Xfc~S„(A) 

k=0 

A 

where S„(A) are Erlang (gamma) laws with pdfs and ch/’s 

fn{x) = Pk{u) = (-r-—'] n = 0 , 1 , 2 ,... 

(n — 1)! “ \\ — lu J 


n = 0,1, 2,... 
n = 0,1, 2 ,... 




N Cufaro Petroni and P Sabino: Cointegrating jumps for energy facilities 


6 


where it is understood that €o = <5o- We will finally denote with N{t) ~ ^(At) and 
M(t) ~ the dependent Poisson processes associated respectively to Tn and 

Sn, and for our purposes we are interested in hnding an explicit form of 


= m, N{t) = n} 

To this end we hrst introduce the parameter 

a/i 


n,m = 0,1 , 2, 


7 = 


A 


and the shorthand notations 


-Kk{a) = e k = 0,l,... 

(3t{n) = -aY £ <n = 0,l,... 


t > 0 


respectively for the distributions of a Poisson ^(a), a binomial 53(n, 1 — a) (it is 
understood that /5o(0) = 1) and a binomial mixture of shifted Erlang laws (^^(A), 
and then we prove (see Cufaro Petroni and Sabino [5]) the following result 


Proposition 2.1. //y > 1, namely ap, > A, we have 
Pm,nit^ ~ 

rn ffL / ■\ / 4 .\ 

Qm,n(s,t) = ^(-1)^ (fc) ^7- V + 1; j + £ + 1; At) 

k=n ^ \ / V / p_n 


0 

Qn,n (^) 

Qm,n (^) Qm,n+l(£^ 


n > m > 0 
m = n > 0 
m > n > 0 


j=k 


e=o 


When 7 < 1, namely ap < X, we have 


^m,n+l(^) T -®m,n(^) 1(^) ^ 7 Tn ^ 0 

Pm,niY) ^ X^n,n{Y) ^n,n+l(^) T Bn,n{Y) T ^n,n{Y) ^ n P 0 

Am,n{t) — ,n+l{t) T Cm,n{Y) Cm, n+lit) m>n>b 


where we define for every n,m > 0 


Am,n(t) = 7lm{pt)'^/3k{n) 

k=0 


1 + 7ik{Xt — apt) — '^^7ij{Xt — apt) 
j=o 


while for n > m > 0, and Xt — apt = w for short, it is 

- — n+l , 

w kl 


Bm,nit) T^mi^pt) ^ ^ J ^ ^ filip T 1) 


k=0 


1=0 


{kp£)\ 


<h £,fc + £ + l, 


1 — a 


-w 
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and for m > n > 1 it is (for n = 0 we have Cm,o(t) = b) 




i=l 


k=n j=0 

[-lY~^~^7rj{Xt)7rra+e-j{apt)^{k + i - j,m + i - j + 1, a/it) 


k + i — j — 1 
k 


and <f>(j + 1; j + £ + 1; At) for 0<£<n<j<m are the confluent hypergeometric 
functions that are in fact elementary functions as proved Cufaro Petroni and Sabino 

H- 

Proof: See Cufaro Petroni and Sabino [5] including the more general case for 
Pm,n{.s,t) = ¥{M{s) = m, N{t) = n}n,m = 0,1,2,... t > 0. Remark that 
in the boundary case 7 = 1 the previous two expressions coherently return the 
same result. □ 


3 The Market Models 

In this section we adapt the model described in Section [2]to the hnancial context. We 
consider an usual Black-Scholes (BS) market and a market with geometric Ornstein 
Uhlenbeck (GOU) processes with jumps similar to the one adopted by Cartea and 
Figueroa [1]. Finally we focus on the Schwartz-Smith model with double jumps that 
can be seen as a complete cointegrated model with jumps. 

Hereafter, compared to Section [21 Ni(t) and N 2 (t) replace M(t) and N(t) and 
Ai and A 2 replace fv and A, respectively. 


3.1 The GBM plus Jumps Case 

Consider a BS market with two risky underlying assets whose dynamics are driven 
by SDEs with the following solution (Merton model): 


Si{T) = exp 


/ 1 A 

log Ri(0) + ( Pi - -an T + aiWi{T) + ^ log 

71 ?.— 1 


, * = 1,2, (3) 


with dWi{t)dW 2 {t) = p^^flt and log-normal jumps: 

,,2 


Ji = Mi exp ( —^ UiZi ) , * = 1,2. 


(4) 


where Zi ~ iV(0,1) and Corr(Z 1 Z 2 ) = We assume that the compound Poisson 
processes and BM are independent. 
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We now concentrate on the logarithm: 

log^,(T) = \ogS,{0)+U^-l^"]T + aiW,{T) + N,{T)\ogM, 


Ni{T) 

+ * = 1 , 2 . 

n,; = l 


(5) 


The equations above can be rewritten as: 

d . ^ / 1 


log^,(T) = log^,(0)+(/r,--a/)T + W(T)logM, 


K 


-fNi(T) + ^/afT + N.{T)i,fH„ i=l,2, 


0 


where given n and m, {Hi, H 2 ) ~ N 

p{j),n,m found in the Appendix lAl 

For simplicity we denote 


1 , P 


(J),n,m 


0 / ’ I ^ 


( 6 ) 

The calculation of 


,y'">(r) = = a^T + n,,f = vf\T) + v. 


{D,n) 


(7) 


where and denote the terminal variances of the continuous and discontin¬ 
uous parts. In case the continuous part of the SDE has a time-dependent volatil¬ 
ity function, it is easy to see that the formulas still hold by replacing vf^\T) by 
f^af(s)ds. 

No-arbitrage conditions imply (see Joshi [3] pag 344): 


fii - r =-AiElJi - 1] i = l,2. 

3.2 The Ornstein-Uhlenbeck plus Jumps Case 


( 8 ) 


Energy markets often display mean-reversion and jumps. We here consider a one- 
factor model plus jumps similar to the one introduced in Cartea and Figueroa [1]. 
Consider a market driven by a stochastic process whose solution is: 

Si(t) = Fi(0,t)exp{l/i(t) + h(t)} , i = l,2, (9) 

where h(t) is a pure deterministic function and Ui{t) is 

.t m) 

Ui{t) = U,{0)e-^^^ + ai V F”* = f/f (t) + f/f (t) (10) 

7o n,=l 

whose SDE is: 


dU,{t) = -kiUi{t)dt + adWi{t) + e-^^%dNi{t). 


( 11 ) 
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yj”* are copies of Yj ~ and Corr{Yi,Y 2 ) = Remark that compared 

to the GBM case the rv's Y[^" are not in terms of logarithms and that our model 
assumes that the size of the jumps decreases in time. The spot SDE is slightly 
different from the one adopted in Cartea and Figueroa [1], indeed the exponential 
term that multiplies the jump component is chosen such that the solution has no 
random jumps with time-dependent jump size. Should we have considered as Cartea 
and Figueroa |T], 

dUiit) = —kiUi{t)dt + aidWiit) + YidNiit), 


the solution would have been 


Ni{t) 


u,{t) = + ai e-’^^^^-^^dWi{s) + ^ Y^'e 


kiT 


ni=l 


where T”* are the jump times. This setting leads to less tractable option formulas 
as it will be shown here below. 

In order to get no-arbitrage conditions, we impose E = F(t, T) and for 


simplicity we look at E[S'(T 
We then need to compute E 
It is well known that: 


= F(0,T) to adjust our parameters and functions. 


uf(t)+uf(t) 


e * 


= E 


P^(t) 


e * 


E 




e * 


E 


Yfd) 


with: 


= exp ( E [Gf (t)] - -Var [U^(t)] ] = 


E[r/f(()] = 

Var[C/f(«)] = R (1 - . 


( 12 ) 


(13) 


Hereafter we will assume that 17*(0) = 0,i = 1,2 that does not change the applica¬ 
bility of the model. Finally we need to calculate: 


E 




/ 


= E 

exp 1 


Ni(t) 

^-kit ^ Y^‘ 

ni=l 




(14) 


Knowing the moment-generating function of the compound Poisson process: 


(j){u) = E 


Ni{t) 


exp < u 




n=l 


= exp{Aii - 1)} 


(15) 


where 


1 
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we easily obtain the required expected value. 


E 




= 0 . 

bi{t) = Xit 


Based on the results above, no-arbitrage is given by hi(t) = —ai{t) — hi(t). 
equations for the spot dynamics above can be rewritten as: 

logS.(() i \ogFi(0,t)-b,{t) + Ni(t)M.e-'“‘ + ^pfe-^'«’N.(t) 


-i (var [t/fti)] ++ 


Var 


uPit) 


Where Hi,i = 1,2 have been dehned in the previous section. 


(17) 

(18) 
The 


(19) 


3.3 The Schwartz-Smith plus Jumps Case 

Consider the two factor Schwartz-Smith model (see Schwartz Smith |13]h 

ft 

Ui{t) = t/i(0)e-^* + cTi / V 

'fo ni=l 

N2{t) 

U2{t) = ^ 2 ( 0 ) fit CT2hh2(f) + ^ ^ 

712 = 1 

U{t) = Ui{t) + U2{t). (20) 

where S{t) = and we assume that the jumps of both process share 

the same distribution Ti,f0 ~ N{M,i>). Simply taking the differential and some 
algebra: 


dU{t) = -k{ia + U2it)-U{t))dt + adW+ Y {e-^^dNi + dN2) (21) 

Should we consider the OU plus compound Poisson as in Cartea and Figueroa P 
the process U{t) would be: 

dU{t) = -k{ia + U2{t)-U{t))dt + aidWi + a2dW2 + Y {dNi + dN2) 

= -k{ia + U2{t)-U{t))dt + adW+ YdN{t). (22) 

where cr^ = af + + 2aia2P^^^- With the latter equation above ([22]), the log of 

the spot process can be expressed in terms of one BM and a compound Poisson-like 
process. 
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With the same procedure outlined in the previous subsection, no arbitrage con¬ 
ditions can be obtained by taking the (conditional) expectation of the spot process: 


E 




= E 


,!7i(0)e-'=*+f7ie-'=(*-'’)dTVi(s)+C/2(0)+/ii+<T2lV2(t)' 


= exp (e [U^{t)] - ^Var (23) 


Assuming once more 17i(0) = 0 and 1 / 2 ( 0 ) = 0 we have: 
E[U^{t)] = pt, 


O’! o . ‘2^pO'i<72 




(1 - e-“) . 


(24) 


For the discontinuous component we have: 


e^’O) = E 




/ 


= E 

exp I 


Ni(t) 


N2it) 


—kt 


E + E 

^1=1 ri.i=l 

ni m 2 

exp ( e-“ 51; n”' + E YT 


+CXD 

^ ^ Pmi,m2 

mi,m2=0 L \ ni=l 

+00 

Pmum2<pY X 0y (1) 

mi,m 2=0 


ni=l 


m 2 


(25) 


As done in Section [5^ no-arbitrage is given by h{t) = —a{t) — h{t) where h{t) can 
be computed numerically. 

In contrast, in the case of Cartea and Figueroa P, we have: 


e^O) = E 




/ 


= E 

exp j 


Ni{t) N2{t) 

ni=l ni=l 


(26) 


that is more complex to treat 

After some algebra, the log-spot dynamics above can be rewritten as: 


logS'j(t) = logF(0, t) — 6(t)-|-/it-|-A'i(t)e [M + e 


^—kt 


V 


PN2{t) M + 


V 


\ [War [U^^\t)] + (e-2^^*Ari(t) + N^it))} + 


'JVar [U^{t)] + (e + N 2 {t)) 


(27) 
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4 Risk-neutral Pricing Formulas 

4.1 The simple European Plain Vanilla Options Case 

In order to simplify the calculation and the notation, we represent the price of a call 
option at time zero c(0) in terms of an abstract BS formula 


c{0) = BS{Po,K,r,T,v,q). (28) 

Where Pq, K, r, T, v, q denote the arguments: initial price, strike, risk-free rate, 
maturity, terminal variance and dividend yield, respectively, for the Black-Scholes 
formula. Following the procedure in Joshi [3] for the Merton model with jumps, 
getting a vanilla option pricing formula will then be a matter of plugging the terminal 
variance, the initial price and the dividend yield into the the abstract formula after 
applying a conditioning argument on the Poisson probabilities. 


• GBM Case. We need to rearrange Equation (jH]) given the market of Equation 
(E]) such that we can apply the abstract BS formula for the GBM-plus-jumps 
case: 


log S,{T) = log ^,(0) + W(T) log M, + A,(l - M,)T - 

The price of a call (put) option on underlying asset i 
jumps market of Equation ([3]) is then: 






+ \ G- 




Hi 


(29) 

= 1 given the GBM-plus- 


c(0) = J2^ni{XiT)BS{S[^\0),K,r,T,v['^’^\T),0). (30) 

n=0 


where 

5{")(0) = 5i(0)Mf exp [AiT(l - Mi)]. (31) 

and is dehned in Equation ([7]). Please note that the formula here is 

slightly different from the one in Joshi |3] pag. 346, because of our dehnition 
of the argument of the abstract BS formula. 

• GOU Case. In contrast for the OU-pIus-jumps market of Equation ([2]), the 
starting point argument for the abstract BS formula is: 

5{")(0) = Fi(0,T)eP^(^\ (32) 


where r;f’"^(T) = Var \uP{T) 


-1- ne 


-2kiT^2 


and 


Piit) 


—bi{t) + ne 


Ip-fcii 

2 




+ Mi 



(33) 
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• Schwartz-Smith Case. Assuming Equation (|20|) a semi-closed form formula 
can be found following the procedure outlined in the GBM and GOU cases. 

OO 

c(0) = =ni,N 2 {T) = n 2 ) 

ni,n 2=0 

(34) 

where 

^(ni,n2)^g^ _ V)(^+V) ^35^ 

and 

[G^(T)] + + TZs) (36) 

Moreover, for both settings in Equations (|20|) and (12^ given the no-arbitrage 
conditions, the simple algorithm to generate the dependent Poisson processes 
of Section |2] offers an easy Monte Garlo implementation to compute the price 
of vanilla options. 

Having obtained the risk-neutral conditions on each underlying asset, it is straight¬ 
forward to obtain formulas for spread options. 

4.2 Spread Options Case 

The application to spread options is the native framework to compare our approach 
with cointegrated jumps compared to other jump-diffusion cases. We start consid¬ 
ering a spread option with zero-strike, based on the results of Margrabe |6], with 
the same conditioning approach applied in Subsection fl4.2p we get that the price of 
a spread option with zero strike given the market of Equation (|3]) or Equation ([9]) 
is: 

OO 

s(0) = = ni;N2{T) = n2) 

ni,n2=0 

(37) 

where ^ ^ 772 ) spread 

terminal variance. The dehnition of S'™ and follows from the subsection above. 

In the literature different analytical approximations are available when the strike 
is not zero (see for instance Deng and Lee |T0] or Kirk |1]), the extension to the jump 
diffusion case is just a matter of adapting the parameters of the approximation. 
Employing Monte Garlo methods is not complicated because the simulation of the 
2-dimensional path is not a complex task and as well as is the 2-dimensional Poisson 
generation. 

^when the Margabe option is written on the spot it can be seen that the price is independent 
of r, in contrast to the situation of a spread option on the forward. 
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We will deserve future studies to a market of more than two assets that have 
a jump component. The current framework cannot cope with multi-asset spread 
options unless one considers that the third asset has no jump term. Pricing multi¬ 
assets spread options then can be tackle via simulation, analytical approximations 
as done in Deng and Lee m and Pellegrino and Sabino im or by applying mo¬ 
ment matching and using one of the solutions available for two legs as explained in 
Pellegrino and Sabino [T6] . 

In the following, we compare three different Poisson models: 

• Independent Jumps. Ni{t) and N 2 {t) are independent Poisson processes. 
The spread option formula is: 

OO 

ni,712=0 

(38) 

• One Common Jump. Ni(t) = N{t) + , i = 1,2, where N(t) and are 

all mutually independent Poisson processes. The spread option formula is: 

CX) 

■5(0) = 7r„i-n(AfT)7r„2-n(A^T)7r„(AT) x 

n=0,ni,n2>n 

• Cointegrated Jumps. Ni{t),i = 1,2 described in Subsection 14.21 The 
spread option formula is: 

CXD 

s(0) = ^{Ni{T) = ni, N2{T) = n2) BS{S['"^\0), 

711,712=0 

(40) 

where Pni,n 2 = P N 2 (T) = ^ 2 ) are dehned in Proposition 12.11 

The payoff of the spread options above considers the values of the two underlying 
at the same time T. Other types of spread options instead look at the two under¬ 
lying at different times, e.g. the payoff may be (S'i(Ti) — S' 2 (T 2 ))’*’, T 2 < Ti. In 
this case one needs to readapt the formulas and consider the probabilities Pni 7 i 2 = 
P(iVi(^i) = N 2 {T 2 ) = 77 , 2 ) and they can be found in Cufaro Petroni and Sabino 

0 . 

5 Numerical Experiments 

In this section we presents the numerical experiments assuming the GBM and GOU 
dynamics plus jumps explained in the previous sections. The case with GBM consid¬ 
ers realistic parameters and is meant to study the spread option values with different 
types of bivariate Poisson processes; in contrast the GOU case is based on real data 
of EEX and Powernext day-ahead prices. 
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Table 1: Parameters of the GBM and Compound Poisson processes 
(a) Continuous Part. (b) Discontinuous Part 



No Jump 

With Jumps 


Case A 

Case B 


5i(0) 

100 

100 

100 

52(0) 

100 

100 

100 

0-1 

0.49 

0.37 

0.2 

0-2 

0.35 

0.23 

0.15 


96 

60 

80 



Case A 

Case B 


99 

50 

Ai 

20 

40 

A 2 

20 

20 

1^1 

0.10 

0.05 

1^2 

0.07 

0.04 

Ml 

1.1 

1.05 

M 2 

1.1 

1.05 


5.1 GBM. Application to Spread Options 

We compare the spread option value obtained using Equations fl38l) - flT0l) changing 
the correlation between the two Poisson processes. In particular, assuming W = 
N{t) + N^{t) we have Cov[Ni{t), N 2 {t)] = Var[iV(f)] = At, then the instantaneous 
correlation is PniN 2 = time-independent; in case of cointegrated jumps 

this can be obtained numerically. Naturally, in case Ai ^ A 2 the perfect correlation 
cannot be obtained. 

We consider two cases: 

Case A . Ai = A 2 = 20, where for cointegrated jumps 7 = a < 1. 

Case B . Ai = 40, A 2 = 20 where for cointegrated jumps 7 < 1 can assume values 
lower or higher than 1, when a < 0.5 and a > 0.5, respectively. 

The parameters are shown in Table [T] In both cases we consider an at-the-money 
spread option with zero strike, K = 0 and maturity T = 1 such that we can use 
the exact Margrabe formula. One can use some approximations techniques for the 
spread option value for non-zero strikes without changing the validity of our tests. 

We also compute, see Table EJ the spread value with a pure GBM with no jumps 
with the parameters of the hrst column of the hrst table in Table [H that are chosen 
such that they match the average spread terminal variance where [•] 

denotes the integer part. 

We hrst observe that the behavior of Equations (l38ll - fj4nil depends on the values 
of the probabilities Pni,n 2 and the values of the BS formulas separately. The former 
quantities do not depend on the distribution of the jumps while the latter ones are 
independent on the structure of dependence between the Poisson processes. 

Figure [1] clearly shows that the expected jump size has a relevant impact in 
the option value because under the assumption Mi = M 2 = 1 the price in almost 
independent on the choice of the Poisson model. 

Figures [2] show the difference among the joint probabilities of the Poisson pro¬ 
cesses. The isolines of the contour plot of Pnin 2 resemble to a sort of ellipse whose 
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(a) Case A (b) Case B 

Figure 1: Spread Option Values in the Cases A and B when Mi = M 2 = 1 

Table 2: Spread Option Values without jumps and independent Compound Poisson 
processes. _ 



No Jump 

Independent Jump 


Case A 

Case B 

Case A 

Case B 

Option Value 

7.27 

11.92 

25.23 

19.27 


axis are parallel to the X-Y axis that is expected for independent Poisson both in 
case A and B. The positive correlation of the Poisson processes in the other con- 
hguration is reflected by the fact that the axis are now rotated counterclockwise. 
In addition for cointegrated Poisson the higher value of the probabilities is more 
concentrated around the expected value. 

This is more evident in the hgure with Ai = A 2 and it is worth noticing once 
more that for 7 > 1 the matrix Pnin 2 is lower triangular for the cointegrated Poisson 
while it is full for the structure with one common Poisson process. Finally, Tabled 
displays the results with expected jump size different from zero. A where chosen such 
that the correlations of the ‘cointegrated’ and ‘common’ Poisson processes coincides. 

The effect of the correlation between the Poisson processes is noticeable: the 
value of the spread option is decreasing when PN 1 N 2 is increasing that is in line with 
the intuition because the spread terminal variance decreases. 

In the case A, the jump sizes are perfectly correlated and the spread option 
values using the common Poisson setting is always higher than the values obtained 
with our methodology. This is somehow reflected by the concentration of the isoline 
of the probabilities. Using a common Poisson reduces the spectrum of jump events, 
for instance in an extreme setting where A = Ai = A 2 , Ni{t) cannot jump more that 
N 2 {t), while this is not the case for the cointegrated Poisson process. In contrast, 
the choice of the Poisson model has no remarkable effect on the price of the spread 
option in the conhguration B. In the parameters conhguration of Tabled! for values 
obtained with the same PniN 2 -i fhe price of the spread option seems to highly depend 
on the number of the jumps of both processes rather than when they occurred and 
that explains the small differences. Furthermore, assuming for instance Af' = 0 
implies A 2 = A and N 2 {t) cannot have more jumps than N 2 {t) that coincides with 
the properties of our model only if 7 > 1, that means that our model gives a richer 
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40 

30 

20 

10 


Independent Poisson 



401— 
30 - 
20 - 
10 - 


10 


15 


Common P§Rson, rlio=0.3’?^ 


30 



35 


40 


40 1 — 
30 - 
20 - 
10 - 


10 


15 


Cointegrated^Poisson, p = 0.§ 


30 


35 


40 


40 1 — 
30 - 
20 - 
10 - 


10 


^Common Pois^Si, p = 0.76, 7 ^ 1 



35 


40 


40 1 — 
30 - 
20 - 
10 - 


10 


d^ntegrated Poi^Pon, rho=0.76,^< 1 


35 


40 


10 


15 


20 


25 


30 


35 


40 


(a) Ai = A 2 = 20 

Independent Poisson 



(b) Ai = 40, A 2 = 20 


r j 1 


1 • j • 


r j 1 
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Table 3: Spread Option Values with Common and Cointegrated Compound Poisson. 



Case A 

Option Value Case A 

Case B 

Option Value Case B 

a 

PNiN2{yo) 

A 

Common 

Cointegrated 

PNiN2{%) 

A 

Common 

Cointegrated 

0.1 

9 

1.80 

24.30 

24.22 

7 

2.09 

18.87 

18.87 

0.15 

14 

2.71 

23.76 

23.64 

11 

3.13 

18.66 

18.67 

0.2 

18 

3.63 

23.20 

23.05 

15 

4.16 

18.45 

18.46 

0.25 

23 

4.55 

22.63 

22.44 

18 

5.19 

18.25 

18.26 

0.3 

27 

5.47 

22.04 

21.81 

22 

6.21 

18.04 

18.05 

0.35 

32 

6.40 

21.42 

21.16 

26 

7.23 

17.83 

17.83 

0.4 

37 

7.34 

20.78 

20.48 

29 

8.24 

17.61 

17.62 

0.45 

41 

8.29 

20.11 

19.78 

33 

9.25 

17.40 

17.40 

0.5 

46 

9.24 

19.41 

19.05 

36 

10.25 

17.18 

17.18 

0.55 

51 

10.20 

18.68 

18.29 

40 

11.25 

16.97 

16.96 

0.6 

56 

11.17 

17.90 

17.49 

43 

12.24 

16.75 

16.74 

0.65 

61 

12.16 

17.08 

16.64 

47 

13.23 

16.53 

16.51 

0.7 

66 

13.15 

16.20 

15.74 

50 

14.21 

16.30 

16.29 

0.75 

71 

14.17 

15.25 

14.78 

54 

15.19 

16.08 

16.06 

0.8 

76 

15.20 

14.21 

13.75 

57 

16.16 

15.85 

15.83 

0.85 

81 

16.26 

13.06 

12.61 

61 

17.13 

15.62 

15.60 

0.9 

87 

17.36 

11.74 

11.33 

64 

18.09 

15.38 

15.37 

0.95 

93 

18.53 

10.14 

9.82 

67 

19.05 

15.15 

15.14 


set of combinations. 

5.2 GOU. Application to Power Interconnectors 

In this section we apply our methodology in order to price a power interconnectors 
between EEX and Powernext modeled as a spread option (Powernext minus EEX). 
We assume that each dynamics behaves as a GOU plus a compound Poisson and cal¬ 
ibrate the parameters. As done for the GBM example we do not consider transport 
costs (no strike) and adopt the Margrabe formula in Equations fl38l) - fHnil . 

In particular, the calculation date is end of December 2015 with a historical 
time window of 2 years for the estimation period. We concentrate then on the 
power spread value for the hrst and second quarters, Ql, Q2 in 201fl§. The choice 
of the maturities above is justihed by the fact that in practice pure spot models are 
used for short maturities, while it is a common practice to use forward models for 
longer maturities and add a volatility premium to the forward volatility. 

Because the expectations in the pricing formulas are in the risk neutral measure, 
some observations are necessary for the estimation of the parameters of the spot 
model. Some approaches propose to assume to live already in the risk-neutral world 
or to incorporate a market price of risk in the drift (see Lucia and Schwartz 0 )- The 

^The technique here discussed does not reflect UGC view. 
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main goal of our work is not the parameters estimation and for sake of simplicity 
we do not consider the market price of risk. 

The estimation procedure can be split into two steps. As a hrst step, after 
hltering out the time-dependent components of each process, one can estimate the 
parameters of the one-dimensional processes, 9i = Xi, As a second 

step then, one can estimate the remaining joint parameters dehned by the two- 
dimensional models 

Consider an equally spaced time grid ... ,tT with tj+i — ti = At and the 
Euler scheme of each SDE in Equation fITT]) 

Ui(t -|- 1) = (1 — kiAt) Ui(t) -|- aiVAtei^t+i + e + l)Tj. (41) 


where 


li{t -\- 1 ) 


1, with probability AjAt 
0, with probability 1 — XiAt. 


(42) 


Hence the transition density is a combination of Gaussian densities: 

Pi [Ui{t + l),i + ^ (1 - + XiAtN' . 

(43) 

Af{x) denotes the density function of a Gaussian random variable and pp(t) = 
(1 — kiAt) Ui{t), niit) = (1 — kiAt) Uiit) + Mie~^'^ and cfit) = ai^/At^ {afit))"^ = 
ufAt + The parameters 6i = (fcj, cTj, A*, Mj, z/j) can be calibrated by mini¬ 

mizing the log-likelihood function with the usual constrains on the parameters: 


T-l 

9i = argmin log {pi {Ui(t + l),t + l\Ui(t),t)). (44) 

t=o 


The calibration of the parameters for the two-dimensional process depends on the 
model specihcation written in Section [21 


• Independent Jumps. In case of independent Poisson processes the joint 
probability are simply: 


Po,o — (1 ~ AiAf) (1 — X2At ), pi^o — AiAt (1 — X2At) 

Po,i = (1 - AiAt) AsAt, Pm = 1 - po,i - pi,o - Po,o- (45) 

The only two remaining parameters to estimate are and and can be 
obtained by minimizing the log-likelihood of the two dimensional process. The 
transition density is 


p{U{t + l),t + l\U{t),t) = 

PP^it)) po,o + AA {p^\t), P^\t)) po,i 
M (/i'^‘^(t), P^^{t)) pi,o M (/i'^'^(t), P-^\t)) pi,i. (46) 

^An example of how to derive the parameters of the 
GOU process in Cartea and Figueroa [1] can be found at 
http://de.mathworks.com/help/fininst/simulating-electricity-prices-with-mean-reversion-and-jump-difFusion.html, 
Although different, the procedure can be adapted to our case. 
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where 


= {Pl{t),P2it)) 

, E^^(t) 

/i^‘^(t) = (pf(t),/i^(t)) 

, E^^(t) 


, E^^(t) 


(ai)^ ) ’ 

{a(r \ 






{a(y p^-^^a(a( \ 
p^-^^a(ai {aif ) ' 


We do not neglect the o (At^) terms that are necessary to estimate p^^\ 


• Common Jumps. One cannot detect the presence of the common Poisson 
process only looking at each log process independently. 

After some algebra, the pair (li(f + 1)) is bi-dimensional Bernoulli rv with: 

Pop = 1 ~ (-^ 1 ^ 'P ^2 + '^) Pi,o = A^At 

Po,i = A^At, = AAt = 1 - po,i - Pip - Po,o-(47) 

In contrast to the case above, we neglect o (At^) terms that implies that we are 
neglecting the possibility that the Poisson processes N^it) and N^if) jump 
simultaneously in the unit of time At. The functional form of the transition 
density is the one of Equation with different probability weights. 


• Cointegrated Jumps 

— Case 1 . 7 > 1 . Based on the results of Proposition fl2.ip up to O(A^) 
terms we have: 

Pup = 1 ~ AiAt, pifi = (Ai — A2)At — Ai — A2^ At^ 

Pop = 0 , Pip = A2At. ( 48 ) 

— Case 2. 0 < 7 < 1 . 

Pop = 1 - (Ai + A2(1 - 7)) At, Pop = A2(1 - 7)At 

Pip = (Ai - 7A2)At, pi,i = 7A2At. ( 49 ) 

Here above we neglect o(At^) terms (see Appendix [B] for the proof of this 
last case) 
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Table 4: Market parameters for EEX and Powernext 
(a) Parameters of the Single Underlyings. 


Market 

k 

o-i 

Pi 

r'i 

Ai 

EEX 

42.50 

1.66 

- 0.10 

0.16 

95.32 

Powernext 

41.64 

1.52 

-0.06 

0.38 

56.74 


(b) Common Parameters 


Method 


A 

a 

Independent 

43 

NA 

NA 

Common 

43 

- 1.0 

NA 

Cointegrated 

43 

NA 

0.44 


Table 5: Powernext-EEX Interconnector Prices 



Interconnector Value (EUR) 


Independent 

Common 

Cointegrated 

Jan 

289.51 

289.64 

310.08 

Feb 

287.40 

287.40 

290.81 

Mar 

263.17 

263.17 

263.37 

Q2 

405.05 

405.05 

405.06 


The results of the calibration are shown in Table 01 The expected jump sizes are 
both negative and their correlation ~ 0 ) is very small and we will neglect it 
hereafter. Comparing the values of A and a or 7 , we can conclude that the correlation 
between the two Poisson processes is of the the order of 45%. In the case of Poisson 
processes with common jumps, the correlation is not time dependent, while in the 
cointegration case the correlation depends on time and it is not a. However, a can 
give a reasonable order of magnitude for the correlation. 

Table [5] shows the values of the power-spread in the quarters Q1 and Q2 and the 
detail of the hrst three months with the different dynamics. As shown in Section 
14.21 the pricing formulas are double sums over probabilities and BS formulas. The 
three models differ in the Poisson probability terms, while the values of the BS 
terms are the same, hence the price difference is only attributable to the different 
bi-dimensional Poisson law. 

The results in Table 0] show that the prices of the interconnector are similar for 
the cases with independent and common jumps for all the maturities. However, 
the model with cointegrated jumps returns higher values for all the maturities with 
a larger difference in January and February. The fact that the prices obtained 
with our methodology are larger than the ones with independent jumps may seem 
counterintuitive and wrong if one relies on linear correlation only, indeed one would 
expect the opposite. The cases of common and independent jumps have similar 
results because the correlation between the Poisson processes of the former Poisson 
choice is not high and apparently the different probability weights do not have a 
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remarkable effect on the prices. In contrast, compared to the other approaches, our 
model gives more emphasis to when the jumps occur and not only to how many 
jumps have occurred. 

A certain news to one of the two markets may cause an instantaneous propagation 
of the information and simultaneous jump in both spot prices (if the jumps belong 
to the common source) while our solution has not this restriction. A same shock 
to one of the two markets may have caused both prices to jumps but at different 
times and that explains why the values are higher. Finally, we observe that at longer 
maturities the difference reduces because of the effect of the mean reversion rate. 

Once more we remark how our methodology, that is parameterized by 7 and a, 
the correlation between the exponential rt>’s that construct the Poisson process, im¬ 
plies a structure that goes beyond the linear correlation and its effects are remarkable 
in the examples illustrated above. 


6 Conclusion and Future Studies 

Based on the concept of self-decomposability we have studied the use of the 2- 
dimensional co-dependent Poisson processes proposed in Cufaro Petroni and Sabino 
[9] to model energy derivatives and in general to price spread options. Due to the 
particular relationships among inter arrival times, we can see this dependence as a 
form of coitegration among jumps. 

To put into the context of modeling energy market and facilities, we have shown 
how to combine 2-dimensional compound Poisson processes with Geometric Browian 
Motions and Geometric Ornstein and Uhlenbeck dynamics. In the latter case, we 
have adopted a dynamics for day-ahead prices that allows (semi-)closed formulas for 
plain vanilla options with an easy derivation of risk-neutral conditions. 

Focusing on the pricing of spread options, we have compared the option prices 
using our methodology and different types of Poisson processes. We have shown 
that our methodology can cope with a wide range of possibilities that go beyond the 
pure correlation between marginal Poisson and can answer several questions that 
arise in the hnancial context. 

In our study we have considered power interconnectors but the applicability can 
be extended to other hnancial situations. Straightforward applications are in credit 
and insurance risk where our approach can answer questions regarding the time of 
contagion or time of propagation of certain information. 

In addition, the self-decomposability and subordination technique can be promis¬ 
ing tools to study dependency beyond the Gaussian-Itd world. For instance, in 
Gufaro Petroni and Sabino [9] we have detailed how to obtain dependent Erlang 
(Gamma) rv^s that can be used to create and simulate dependent variance gamma 
processes. Furthermore in recent papers, Sexton and Hanzon [T3] and Hanzon et ah 
[2] have studied the use of two sided Exponential-Polynomial-Trigonometric (ETP) 
density functions to option pricing where EPT are distributions with a strictly proper 
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rational characteristic function. Due to the fact that the Erlang and exponential 
distributions belong to this class, it will be worthwhile to investigate the use of 
self-decomposability to create dependence for this larger class of distributions. 
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A Calculation of 

• The GBM plus Jumps Case. 

piJ),n,m ^ [aifTi(T) + a2W2{T) + = 

= =- / -• (50) 

• The Ornstein-Uhlenbeck plus Jumps Case. 

^{J),n,m _ Corr{Li,L2) = 

= Jblhhl - - -(51) 

yr;f’”^(T) 4 '^’™^(T) 

where Li = ui e~^^^^~'^'>dWi{s)+^/nuie~’^^^Zi and L 2 = a 2 Jq e~’^^^^~^'>dW 2 {s) + 
^/mu2e~’^^^Z2. 


B Calculation of ] 9 o,o, 1^0,1, Pip and pi^i 


Based on the results in Cufaro Petroni and Sabino the joint cdf of Xi ~ Ci(Ai) 
and X 2 (£i(A 2 ) is 


H{xi,X2) 



(52) 


For At small we can assume that no more that one jump can occur hence: 


Pi,i = P(Xi < At, X2 < At) = H{At, At) 

Pi,o = P(Xi < At, X 2 > At) = Fi{At) - H{At, At) 

Po,i = P(^i > At, X 2 < At) = F 2 {At) - if (At, At) 

Po,o = P(Xi > At, X2> At) = 1- pi,i - po,i - Po,i 


Finally, the results in section 15.21 are obtained considering only 0{At) terms and 
splitting between 7 > 1 and 0 < 7 < 1 . 
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